SPARTA, Inc. a Parsons Company 

Lake Forest, CA 92630-8873 


Final Technical Report for FA9550-09-C-0207 

"Partitioning of Electromagnetic Energy 
Inputs to the Thermosphere during 
Geomagnetic Disturbances" 

June 2012 
Prepared by: 

Keith Siebert, Principal Investigator 
Applied Research Associates (ARA) 

39 Simon Street, #15 
Nashua, NH 03060 


Prepared for: 





REPORT DOCUMENTATION PAGE 

Form Approved 

OMBNo. 0704-0188 

data needed^ and completing and reviewing this collection of information. Send comments regarding this burden estimate or any other aspect of this collection of information, including suggestions for reducing 
this burden to Department of Defense, Washington Headquarters Services, Directorate for Information Operations and Reports (0704-0188), 1215 Jefferson Davis Highway, Suite 1204, Arlington, VA 22202- 

vaHd OM^confro^rfumber d PLEASE DO^O^RCTIHRN YOUIRDRM^O TH e!^BOVE /UDDRESS.' ^ ^ ^ ^ ^ C ° mPly ^ * C ° lleCtl ° n ° f informatlon ,f ll does not dlSplay a current| y 

1. REPORT DATE (DD-MM-YYYY) 1 2. REPORT TYPE 

30-06-2012 | Final 

3. DATES COVERED (From - To) 

30-09-2009 - 30-06-2012 

4. TITLE AND SUBTITLE 

Final Technical Report for FA9550-09-C-0207 

Partitioning of electromagnetic energy inputs ihto the 
thermosphere during geomagnetic disturbances 

5a. CONTRACT NUMBER 

FA9550-09-C-0207 

5b. GRANT NUMBER 

5c. PROGRAM ELEMENT NUMBER 

6. AUTHOR(S) 

Siebert, Keith D 

5d. PROJECT NUMBER 

5e. TASK NUMBER 

5f. WORK UNIT NUMBER 

7. PERFORMING ORGANIZATION NAME(S) AND ADDRESS(ES) 

Applied Research Associates 

39 Simon Street, #15 

Nashua, NH 03060 

8. PERFORMING ORGANIZATION REPORT 
NUMBER 

9. SPONSORING / MONITORING AGENCY NAME(S) AND ADDRESS(ES) 

AFOSR/NS (Dr. Kent Miller) 

875 North Randolph Street 

Suite 325 f Room 3112 

Arlington, VA 22203 

10. SPONSOR/MONITOR’S ACRONYM(S) 

11. SPONSOR/MONITOR’S REPORT 
NUMBER(S) 

AFRL-OSR-VA-TR-2012-1160 

12. DISTRIBUTION / AVAILABILITY STATEMENT 

DISTRIBUTION A: APPROVED FOR PUBLIC RELEASE 

13. SUPPLEMENTARY NOTES 

14. ABSTRACT 

We investigate a number of fundamental processes that govern thermospheric response to 
electromagnetic energy inputs. Using a sophisticated numerical model of the coupled solar 
wind - magnetosphere - ionosphere - thermosphere, we focus oij three main subject areas: 
Poynting flux energy partitioning, horizontal neutral density structures, and feedback of 
neutral dynamics to the coupled system. Simulations show that approximately 90 per cent of 
the Poynting flux into the thermosphere goes into Joule heating under steady state 
conditions. Under southward IMF orientations, local height-integrated Joule heating rates 
overestimate PoynKifig flux in regions of high conductance, and underestimate Poyntihg flux in 
regions of lower conductance, confirming the prediction of Richmond (2010). At E-region 
altitudes, neutral density cell structure to is correlated with: divergences of the J X;’^. 
force. These divergences arise because of Hall conductivity effects which rotate current- 
density and ion-velocity vectors from their non-divergent F-region patterns. 

15. SUBJECT TERMS 


16. SECURITY CLASSIFICATION OF: 


a. REPORT 


b. ABSTRACT 


c. THIS PAGE 


17. LIMITATION 
OF ABSTRACT 

SAR 


18. NUMBER 19a. NAME OF RESPONSIBLE PERSON 
OF PAGES Keith D. Siebert 

19 19b. TELEPHONE NUMBER (include area 

(603) 595-2139 


Standard Form 298 (Rev. 8-98) 

Prescribed by ANSI Std. Z39.18 







Table of Contents 


1. Summary.1 

2. Introduction.1 

2.1. Poynting Flux Energy Partitioning.1 

2.2. Investigation Objectives.3 

2.3. Organization of Report.4 

3. Methods, Assumptions and Procedures.4 

3.1. Integrated Space Weather Prediction Model.4 

3.2. Simulation Methods.4 

4. Results and Discussion.6 

4.1. Dependence of Energy Partitioning on IMF Clock Angle and Season.6 

4.2. Effect of Drag-Driven Neutral Winds on Thermospheric Density.10 

4.3. Feedback Mechanisms.15 

5. Conclusions.16 

6. Personnel Supported.17 

7. Publications.18 

8. Interactions/Transitions.18 

9. References.18 

List of Figures 

Figure 1. Pedersen and Hall conductance maps for (a) vertical geographic pole 
(northern hemisphere); (b) solsticial pole (northern hemisphere); (c) 
solsticial pole (southern hemisphere, viewed through the northern 
hemisphere).6 

Figure 2. Joule dissipation rate (left) and mechanical work rate (right) at 120 km 

altitude, northern polar region, IMF South. Red star in left plot indicates 
intersection point of radial line used for 1-D plot in Figure 3.7 

Figure 3. Joule heating rate as a function of altitude in units of energy density per 

time (top) and specific energy per time (bottom).8 























Figure 4. Poynting flux maps (left) and height-integrated power maps (right) in 

northern polar region for (a) South IMF; (b) East IMF; (c) North IMF. On each 
map the "+" symbol locates peak value.9 

Figure 5. Ion velocity vectors at 400 km (left) and 120 km (right) for (a) IMF South; (b) 

IMF East; (c) IMF North. Vectors are color-coded according to magnitude 

with black indicating speeds equal or greater than 1 km/s for (a) and (b), and 

0.5 km/s for (c).11 

Figure 6. Neutral velocity vectors in the rotational frame (top) and neutral density 
contours (bottom) at 400 km. From left to right: IMF East - northern 
hemisphere; IMF East - southern hemisphere; IMF South; IMF North. 

Vectors are color-coded according to magnitude with black indicating speeds 
equal or greater than 1 km/s.12 

Figure 7. Figure 6 from Kwak et al. (2009). Difference densities at 400 km altitude over 

the southern hemisphere for different IMF orientations.13 

Figure 8. Neutral velocity vectors in the rotational frame at 120 km. From left to right: 

IMF East - northern hemisphere; IMF East - southern hemisphere; IMF 
South; IMF North. Vectors are color-coded according to magnitude with 
black indicating speeds equal or greater than 0.1 km/s.14 

Figure 9. J X B vectors (top) and neutral density contours (bottom) at 120 km. From 
left to right: IMF East - northern hemisphere; IMF East - southern 
hemisphere; IMF South; IMF North. Vectors are color-coded according to 
magnitude with black indicating J X B force magnitudes equal or greater 
than 1 X 10 11 dyne/cm 3 .14 

Figure 10. Electrodynamic quantities for simulations with zero neutral velocity (left 

column) and steady-state neutral velocity (right column), (a) J X B vectors at 
300 km altitude; (b) Parallel current density contours at 500 km altitude 
(orange-red indicates current into ionosphere; blue indicates current out of 
the ionosphere); (c) J X B vectors at 120 km altitude.17 

List of Tables 

Table 1. Simulation Parameters.5 

Table 2. Electromagnetic Power Energy Partitioning.8 












Partitioning of Electromagnetic Energy Inputs to the Thermosphere 
Final Report for FA9550-09-C-0207 


1. Summary 


We investigate a number of fundamental processes that govern thermospheric response to 
electromagnetic energy inputs. Using a sophisticated numerical model of the coupled solar wind - 
magnetosphere - ionosphere - thermosphere, we focus on three main subject areas: Poynting flux 
energy partitioning, horizontal neutral density structures, and feedback of neutral dynamics to the 
coupled system. 

Numerical simulations show that approximately 90 per cent of the Poynting flux into the thermosphere 
goes into Joule heating under steady-state conditions. Thus in a global sense, measurements of Joule 
heating rates provide a good approximation to total electromagnetic energy input; however, locally this 
is not always the case. Under southward IMF orientations, local height-integrated Joule heating rates 
overestimate Poynting flux in regions of high conductance, and underestimate Poynting flux in regions 
of lower conductance. This confirms the prediction of Richmond (2010). 

We confirmed the analysis of Schoendorf et al (1996) that attributed the existence of anti-cyclonic low 
density cells to centrifugal effects in the presence of large drag-driven neutral flows. At E-region 
altitudes, we present a new theory that relates neutral density cell structure to divergences of the ] xB 
force. At E-region altitudes these divergences arise because of Hall conductivity effects which rotate 
current-density and ion-velocity vectors from their non-divergent F-region patterns. 

We find that increases in neutral wind speeds under steady driving conditions do not significantly alter 
electromagnetic energy flows into the ionosphere / thermosphere. While F-region closure currents are 
significantly reduced, the much larger E-region closure currents are not, owing to the much larger 
neutral densities, and thus much greater coupling times. Thus parallel current densities and Poynting 
fluxes into the ionosphere / thermosphere are not greatly altered by neutral dynamics. 


2. Introduction 


The dynamics of Earth's thermosphere represents a key challenge in the area of SSA. Aerodynamic drag 
is the largest source of uncertainty in precision orbit determination for satellites operating below about 
600 km altitude. Aerospace missions affected by drag errors include satellite orbit location and 
prediction, collision avoidance warnings, reentry prediction, lifetime estimates, and attitude dynamics 
(Marcos, 2006). Because aerodynamic drag is dependent on thermospheric density, nowcasting and 
forecasting of the state of the thermosphere is necessary for reduction of drag errors. 

2.1. Poynting Flux Energy Partitioning 

The thermosphere is one element of a complex coupled system that includes the ionosphere, 
magnetosphere and the sun. Thermospheric density responds to energy inputs generated within this 
system. The sun is the ultimate driver, and provides energy in the form of EUV radiation and solar wind 
kinetic energy. Solar EUV radiation directly heats and ionizes the thermosphere. Solar wind energy is 
transferred indirectly through an interaction with Earth's magnetic field that gives rise to an 
electromagnetic dynamo. Dynamo energy flowing into the thermosphere is mediated by the 
ionosphere. Unlike solar EUV radiation to which the thermosphere responds passively (i.e. without 
feedback), the thermosphere is a reactive element in the coupled system through which solar wind 
energy flows. Thus the thermospheric response can alter energy flow through the system. It is this 
aspect of thermospheric dynamics that presents the greatest challenge to modeling efforts because it 
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requires a treatment of the solar wind - magnetosphere - ionosphere - thermosphere system in its 
entirety. 

The two primary sources for energy transfer from the solar wind dynamo to the thermosphere are 
particle precipitation and electromagnetic field energy or Poynting flux. This report focuses on the 
response of the thermosphere to Poynting flux. While particle precipitation does produce heating of the 
thermosphere, in general the total power input from Poynting flux is greater than that of particle 
precipitation by roughly a factor of two (Waters, 2004). In our study, particle precipitation is 
investigated primarily for its role in altering ionospheric conductance profiles. Variations in ionospheric 
conductance alter the total power delivered to the thermosphere by Poynting flux. 

The mechanism by which Poynting flux is transferred to the thermosphere is most easily understood by 
considering the thermosphere and ionosphere to be separate hydrodynamic and magnetohydrodynamic 
(MHD) fluids respectively. Convection of the ionosphere fluid associated with the Poynting flux electric 
field exerts a drag force on the thermosphere fluid to which it is collisionally coupled. This accelerates 
the thermosphere and heats it frictionally. In addition to momentum transfer, collisional processes 
transfer heat between the ion and neutral fluids. 

Polar Poynting flux distributions depend on IMF orientation. Using DE2 data, Gary et al. (1995) showed 
that for southward IMF Poynting flux was greatest in the auroral zones at dawn and dusk. For northward 
IMF, they found that the Poynting flux was greatest in the polar cap near noon. Waters et al. (2004) also 
showed Poynting flux peaking in the auroral zones for southward IMF, but noted that the Poynting flux 
was not uniform in local time and tended to maximize in "hot spots" whose location and extent differs 
between hemispheres and depends on the IMF. They also found that the total Poynting flux power 
integrated over the polar cap is an appreciable fraction of the total electromagnetic energy, accounting 
for up to half as much as the auroral zone input. 

The paradigm for partitioning Poynting-flux energy into thermal and mechanical energy pools is based 
on Poynting's theorem, which for time-independent fields can be written as (e.g. Thayer, 1995): 

- V ■ S = - V ■ £ xSB = J-E = J-W + U -J xB 

Flere S is the Poynting flux, SB is the perturbation magnetic field (deviation from Earth's background 
dipole field), / is the electric current density, E the electric field in the inertial frame, E' is the electric in 
the frame of the moving neutral fluid, and U is the neutral fluid velocity. The / ■ E' term on the far right 
hand side is typically identified as Joule dissipation that goes directly into thermal energy. However, care 
must be taken with this interpretation. Vasyliunas and Song (2005) argue that "Joule dissipation" is a 
misnomer and this heating term is more accurately described as frictional heating arising from the 
relative motion of plasma and neutral fluids. This heating causes the thermosphere to expand upward 
and thus produces density enhancements along satellite orbit tracks. 

The U ■] xB term represents the portion of the electromagnetic energy that goes into mechanical 
work done on the neutral wind. When the ionosphere is near equilibrium a force balance exists between 
the neutral-ion drag force (i.e. the drag force exerted by the neutrals on the ions) and the ] xB force: 

PiVinOK - i£) = / x 5 

Here pfis the ion mass density, v in is the ion-neutral collision frequency, and v[ and v £ are the ion and 
neutral velocities respectively. By momentum conservation, the ion-neutral drag force is exactly equal to 
the neutral-ion drag force (i.e. the drag force exerted by the ions on the neutrals). In this way the 
] xB force is directly transferred to the neutrals. The mechanical work term can thus be viewed as the 
work done by the ] xB force on the neutrals. Unlike the Joule heating term, the mechanical work term 
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can be positive or negative. Positive work indicates energy transfer from ions to neutrals; negative work 
indicates energy transfer from neutrals to ions. 

The presence of the neutral wind adds significant complexity to the energy partitioning problem, and is 
the reason it has become a subject of detailed study (Thayer et al., 1995, 1998, 2000). These studies 
have shown that ignoring it can result in significant errors in the determination of global thermospheric 
heating rates. Furthermore, for given solar wind conditions, electromagnetic energy partitioning 
depends on spatial location. Over the polar cap Poynting flux is more heavily weighted toward the 
mechanical energy pool, while auroral-zone Poynting flux is weighted toward thermal (Waters et al., 
2004). 

Like Joule dissipation, electromagnetic energy that goes into the mechanical energy of neutral winds can 
also alter thermospheric density profiles. The ultimate response of the thermosphere is governed by the 
interplay of the drag force with other neutral forces such as pressure gradients, centrifugal and Coriolis 
forces, and viscosity. Illustrating this point, Schoendorf et al. (1996) describe a mechanism by which the 
neutral winds drive the high-latitude density structure. 

The mechanism is initiated by ion drag forcing high-latitude neutral winds into a cellular circulation 
pattern similar to the ion velocity pattern. They find that for slower neutral circulation speeds, the 
thermosphere behaves as expected for typical geostrophic flow: pressure gradient and Coriolis forces 
oppose each other. Therefore, in the dusk sector where flow is typically anticyclonic, neutral circulation 
encompasses a high density region while in the dawn sector where flow is typically cyclonic, neutral 
circulation encompasses a low density region. When situations occur in which ion drag causes neutral 
velocities to be large, thermospheric density behavior on the dawn side remains the same, but in the 
anticyclonic dusk sector density is depleted instead of enhanced. The mechanism for this depletion was 
posited to be an anomalous "antibaric" low for which the pressure gradient and Coriolis forces act in the 
same direction, towards the center of the circulation cell. Schoendorf et al. also determined that high 
latitude Joule heating is responsible for the general increase in neutral densities as geomagnetic 
conditions strengthen, while the aforementioned neutral circulation reduces densities around 
circulation centers. The combination of these influences produces the structure seen in high-latitude 
neutral mass densities. Using CHAMP data, Kwak et al. (2009) came to a similar conclusion when they 
investigated high-latitude neutral densities for various IMF conditions, and determined that ionospheric 
convection drives the thermospheric winds which influence the high latitude neutral densities, and also 
that auroral and cusp region density variations are influenced by local heating associated with 
ionospheric currents. 

The thermosphere is an active element of the coupled system. Electromagnetic energy input into the 
thermosphere changes properties of the thermosphere that alter subsequent energy partitioning, and 
can alter the Poynting flux itself as the global electric circuit reconfigures to accommodate changes to its 
fundamental elements. Such changes are generally on longer time scales, and are complicated by the 
different response times of electric field, conductivity, and neutral wind. Determination of the operative 
feedback loops within the coupled system is critical for understanding its evolution and for developing 
accurate predictive models. 

2.2. Investigation Objectives 

Our investigation was organized by four different topic areas: 

1) Magnetic field effects: How does IMF orientation and dipole tilt affect electromagnetic energy input 
and partitioning into the thermosphere? 
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2) Ion drag effects: How does the rotation of the ion velocity pattern with altitude (due to Hall 
conductivity) affect neutral density distributions? 

3) Feedback mechanisms: How do thermospheric dynamics alter electromagnetic energy input and 
partitioning? 

2.3. Organization of Report 

Section 3 contains a brief description of the simulation model, and the input parameters that specified 
the set of simulations we performed for this investigation. We also describe diagnostics we used to 
analyze simulation data sets. Section 4 is divided into three sections, with each section devoted to one 
of the topic areas discussed above. Section 5 contains the main conclusions from this investigation. 


3.1. Integrated Space Weather Prediction Model 

The simulation tool used in this investigation was the Integrated Space Weather Prediction Model (ISM). 
ISM is a global model of the coupled ionosphere - thermosphere - magnetosphere system. It has been 
developed under sponsorship of the Defense Threat Reduction Agency (DTRA) by a technical team that 
comprises industry (SPARTA and ARA) and academia (Boston University, Dartmouth College, and Rice 
University). AFRL has also participated during ISM development in an advisory capacity. ISM has been 
used in a variety of scientific investigations sponsored by NASA and NSF. 

ISM is based on the equations of Two-Fluid MHD (White, et al., 2001). The two fluids are a plasma fluid 
(representing the ion-electron plasma of the ionosphere, magnetosphere, and solar wind) and a neutral 
fluid (representing the atomic and molecular species of the thermosphere). The plasma and neutral 
fluids are coupled collisionally (momentum and energy exchange) and chemically. In regions where the 
neutral density goes to zero, Two-Fluid MHD reduces to the standard MHD equations. 

ISM advances the Two-Fluid MHD equations on a grid system that extends continuously from the 
"bottom" of the ionosphere (typically set at 80 km) upward through the magnetosphere and into the 
solar wind. The outer grid boundaries are typically set at 30 R E sunward, 300 R E tailward, and 50 R E dawn 
and duskward. Since the MHD equations are essentially statements of the conservation of mass, 
momentum and energy, ISM's use of a computational grid that continuously spans the relevant physical 
volume, combined with a numerical scheme that preserves the conservation properties of the 
equations, makes it well-suited for studying energy flow within the solar wind - magnetosphere - 
ionosphere - thermosphere system. 

Conductivities are not explicitly computed in ISM except as outputs for diagnostic purposes. Rather, the 
ion-neutral coupling is expressed in terms of collision frequencies (momentum coupling) and 
heat-transfer coefficients (energy coupling). The approach is consistent with the formulation of 
Vasyliunas and Song (2005) and permits expression of Joule heating in terms of the frictional heating 
term, PiVin\Vi-Vn\ 2 . 


3.2. Simulation Methods 

Parametric Variations 

Table 1 lists parametric variations performed for this study. All simulations had a common solar wind 
velocity and density of 400 km/s and 5 cm 3 respectively. A rotating Earth was used with the 
geomagnetic and geographic poles aligned. Simulations were performed in three phases. First a quasi- 
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steady magnetosphere was modeled with a static thermosphere to establish the ionospheric driver. This 
driver was then used to "spin-up" the thermosphere for a simulation time interval of 24 hours, after 
which the thermosphere was also quasi-steady. Thereafter simulations were fully coupled with fully- 
dynamic magnetosphere, ionosphere, and thermosphere regions. 


Table 1. Simulation Parameters 


Case 

IMF B x 
(nT) 

IMF B y 
(nT) 

IMF B z 
(nT) 

Conductance 

Orientation of 
Geographic Pole 

Moderate Driver 

0 . 

0 . 

-2. 

Uniform 

(E P = Z H = 1 mho) 

Vertical 

IMF South 

0 . 

0 . 

-5. 

Solar + Auroral 

Vertical 

IMF East 

0 . 

5. 

0 . 

Solar + Auroral 

Vertical 

IMF North 

0 . 

0 . 

5. 

Solar + Auroral 

Vertical 

Summer Solstice 

0 . 

0 . 

-5. 

Solar + Auroral 

Northern Hemisphere 
Summer Solstice 

Strong Driver 

0 . 

0 . 

-20. 

Solar + Auroral 

Vertical 

IMF Rotation 

0 . 

0 . 

-5. To +5. 

Solar + Auroral 

Vertical 


Figure 1 shows Hall and Pedersen conductance profiles for the cases labeled "Solar + Auroral" in the 
conductance column of Table 1. These conductance profiles are not specified directly; rather, 
conductances are determined by specification of solar and auroral ionization sources. Figure 1 contains 
separate plots for northern and southern hemispheres for the "Summer Solstice" case, since 
conductances are not symmetric for that case. 

Diagnostics 

We used the ISM_VIEW visualization application to generate contour plots and vector streamline plots 
of both native ISM cell quantities as well as derived quantities such as Poynting flux. Because ISM carries 
the perturbation magnetic field, SB, as a cell quantity throughout its computational volume, we are able 
to generate maps of the Poynting flux, E x SB, on altitude surfaces at the top of the thermosphere. We 
used separate modules to integrate the Poynting flux over this surface to compute the total 
electromagnetic energy input into the thermosphere. We also created modules to compute the volume 
integrals of the Joule heating rate and the mechanical work rate. With these diagnostics we were able to 
compute both the total energy input and energy partitioning for all of the cases listed in Table 1. 
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Pedersen 


Hall 



i (mhos) 


Figure 1. Pedersen and Hall conductance maps for (a) vertical geographic pole (northern hemisphere); 
(b) solsticial pole (northern hemisphere); (c) solsticial pole (southern hemisphere, viewed through the 
northern hemisphere) 


4. Results and Discussion 


4.1. Dependence of Energy Partitioning on IMF Clock Angle and Season 

In this section we look at the spatial distribution of Poynting flux, Joule dissipation, and mechanical work 
rate, as well as height, area, and volume integrals of these quantities over the poles. One of the goals of 
this section is to assess the degree to which Joule dissipation rates can be used in data analyses and 
model analyses as a surrogate for Poynting flux. Figure 2 contains north-pole contour plots of the Joule 
dissipation rate and mechanical work rate at 120 km altitude for the IMF South case. Note that 
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mechanical work rate magnitudes are significantly smaller than Joule dissipation magnitudes. (The color 
scale for the mechanical work rate plot is a factor of 10 smaller than the scale for the Joule dissipation 
plot.) Both quantities are peaked in auroral conductance regions, with lesser maxima in the convection 
throat region. Joule dissipation is a positive definite quantity, while mechanical work rate is a signed 
quantity. In the mechanical work rate plot, tan colors indicate positive values and blue indicates 
negative values. The convection throat region is a region of negative mechanical work rate, indicating 
that neutrals are doing work on the ions in this region; however, the magnitude of this work rate is 
about a factor of 40 smaller than the Joule dissipation rate. 

Figure 3 shows the altitude distribution of the Joule heating rate along a radial line through the red star 
superimposed on the left plot of Figure 2. The top panel is the volumetric Joule heating rate in units of 
erg/cm 3 -s. It exhibits a single peak at approximately 120 km altitude (the altitude of the contour plots in 
Figure 2). This is not surprising since 120 km is an E-region altitude where conductivities are peaked, and 
Joule heating goes like cr|£"| 2 , where a is the conductivity and E' is the electric in the frame of the 
moving neutrals. The bottom panel of Figure 3 is the specific Joule heating rate in units of erg/g-s. It 
corresponds to the time rate of change of neutral temperature, and exhibits both an E-region peak and 
a much larger F-region peak. The F-region is a local maximum of Pedersen conductance, and although 
conductance is smaller than in the E-region, the lower neutral density there results in a larger specific 
heating rate. 


Joule Dissipation Rate 



Mechanica I Work Rate 



Figure 2. Joule dissipation rate (left) and mechanical work rate (right) at 120 km altitude, 
northern polar region, IMF South. Red star in left plot indicates intersection point of 
radial line used for 1-D plot in Figure 3. 


From Poynting's theorem we know that in steady-state the volume integral of the sum of Joule heating 
and mechanical work rates should be equal to the surface integral of Poynting flux normal to an altitude 
surface above the thermosphere. This assumes that we choose a volume whose latitudinal extent away 
from the pole is great enough that horizontal Poynting flux is zero at the latitudinal boundary. It also 
assumes that vertical Poynting flux is zero at the volume's low-altitude boundary. Because our 
simulations allow us to compute Poynting flux independently from Joule heating and mechanical work 
rates, we can assess the degree to which they satisfy Poynting's theorem. Table 2 contains the relevant 
integrated quantities for the IMF South, East, and North simulations. The second column is the Poynting 
flux integrated over a 500 km altitude surface. The values in this column indicate that the total 
electromagnetic energy input into the thermosphere decreases as the IMF is rotated from south to 
north. This is consistent with the findings of other studies (e.g. Weimer, 2004). The third and fourth 
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columns are the volume integrals of Joule heating and mechanical work rates respectively. The fifth 
column is the sum of columns three and four. This number should be approximately equal to the 
integrated Poynting flux in column two. The fractional difference between it and integrated Poynting 
flux is given in column six. The small values here confirm simulation operation. Exact agreement is not 
expected given that the simulation is not in an exact steady state. 



SM LT 10.50 10.50 10.50 10.50 10.50 10.50 10.50 10.50 10.50 

SM AL T 100, 150. 200. 25 0. 300. 350. 400. 450. 500- 

Figure 3. Joule heating rate as a function of altitude in units of energy density 
per time (top) and specific energy per time (bottom). 


Given that simulation results agree with expectations from Poynting's theorem, we can now assess the 
energy partitioning of Poynting flux into thermal and mechanical energy pools. Column 7 of Table 2 
contains the fraction of Poynting flux going into mechanical work. For all IMF orientations shown, 
mechanical work is on the order of 10 per cent of the total electromagnetic energy going into the 
thermosphere. This relatively small value indicates that at least in an integrated sense, Joule heating 
rates account for almost all of the electromagnetic energy input into the ionosphere / thermosphere. 

Table 2. Electromagnetic Power Energy Partitioning 


Case 

Poynting 

Power 

(GW) 

Joule 

Heating 

(GW) 

Mechanical 

Work 

(GW) 

Joule + 

Mechanical 

(GW) 

Fractional 

Difference 

Mechanical 

Work 

Fraction 

IMF South 

67 

57 

7.8 

65 

3% 

12% 

IMF East 

34 

28 

4.3 

32 

4% 

13% 

IMF North 

3.7 

3.3 

0.3 

3.6 

3% 

8% 


Next we want to examine whether Joule heating rates correspond to electromagnetic energy input in a 
local sense. To do this we need to compare Poynting flux maps generated on an altitude surface above 
the thermosphere with maps of field-line integrated Joule heating rates. The degree to which they agree 
represents the degree to which Poynting flux at a particular point above the thermosphere is dissipated 
locally in the flux tube volume below that point, and thus the degree to which field-line integrated Joule 
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heating can approximate local Poynting flux and vice versa. Work by Richmond (2010) showed that this 
equality can only be satisfied if the side boundary of a local flux tube volume is an equipotential. Figure 
4 contains maps of Poynting flux normal to a 500 km altitude surface and maps of height-integrated 
Joule heating rate for the IMF South, East, and North cases. We have chosen height-integration rather 
than field-line integration for computational simplicity. For the polar regions under consideration, where 
field lines are nearly vertical, this should be a good approximation. The color scale for these plots has 
both negative and positive values since we are plotting a single Poynting flux component (the 
component normal to the constant-altitude plotting surface), and it is a signed quantity. The convention 
we use is positive values (tan) indicate flux into the thermosphere, and negative values (blue) indicate 
flux out of the thermosphere. Figure 4 shows that for all of these steady-state cases, essentially all of the 
Poynting flux is directed downward into the thermosphere. 



Figure 4. Poynting flux maps (left) and height-integrated power maps (right) in northern polar region 
for (a) South IMF; (b) East IMF; (c) North IMF. On each map the "+" symbol locates peak value. 
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Figure 4 shows that for the IMF East and IMF North cases, Poynting flux maps are in reasonable 
agreement with maps of height-integrated Joule heating. Peak values are in approximately the same 
location (indicated by the "+" symbol overlaid on the plots) and agree to better than 20 per cent. For the 
IMF South case the maps are in poorer agreement. There, Poynting flux is peaked in the convection 
throat region, while the integrated Joule dissipation rate is peaked in the auroral region. These results 
are consistent with one of the main conclusions of Richmond (2010). Fie found that the Poynting flux 
underestimates electromagnetic energy dissipation in regions of high ionospheric Pedersen 
conductance, and overestimates dissipation in regions of low conductance. The top two plots of Figure 4 
confirm this with the Poynting flux being less than the dissipation rate in the auroral ring (region of high 
conductance) and greater than the dissipation rate in the convection throat (region of low 
conductance). 


4.2. Effect of Drag-Driven Neutral Winds on Thermospheric Density 

In this section we investigate how mechanical work done by the / x B force (or equivalently the ion 
drag force) on the neutrals alters thermospheric density profiles. Ultimately this is an investigation into 
how that force, together with other thermospheric forces such as pressure gradient, centrifugal, and 
Coriolis forces, determines particular steady-state horizontal density structures at different altitudes. 
Flere we focus on two altitudes that capture two distinct ion velocity patterns: 400 km and 120 km. 
Figure 2 illustrates these patterns for the IMF South, East, and North cases. At the higher altitude the 
flow patterns are largely rotational. IMF South exhibits two circulation cells of approximately equal 
magnitude (the tilt of the pattern is due to the Flail conductance). IMF East also has two circulation cells, 
but the poleward cell is considerably larger than the equatorward cell. IMF North has four circulation 
cells with the largest pair having a sunward convection region at its center. These are the canonical 
patterns for their respective IMF orientation (e.g. Weimer, 2004). In contrast, at the lower altitude 
where the Flail conductivity dominates, the patterns are largely irrotational. Note the regions of 
convergence and divergence of the velocity vectors in these patterns. We can expect that the 
compressibility properties of these patterns will be imposed on the neutral fluid through the ion drag 
force. 

Neutral Flows and Density at 400 km 

Figure 3 shows the neutral velocity flow pattern and horizontal density profile of the thermosphere at 
400 km altitude for the IMF East, South, and North cases. Two sets of plots are shown for the IMF East 
case because of the anti-symmetry of the northern and southern polar regions for this case. In all cases 
the neutral flow vectors track pretty closely the corresponding ion velocity vectors (see left column of 
plots of Figure 2). Because the neutrals are not constrained by the magnetic field, the momentum 
imparted to them by the ions at high latitudes enables them to slip through the ions at lower latitudes. 
This is particularly evident in the IMF East case where the low-latitude, nightside vectors are directed 
toward lower latitudes rather than turning sunward as they do in the corresponding ion vector plot. 

The color coding on the neutral density plots of Figure 3 is such that blue indicates a density low, and 
orange-red indicates a density high. In all cases a low density cell is centered on the strongest circulation 
centers. For the IMF East case there is only one low density cell in each hemisphere. For the IMF North 
and South cases there are two low density cells. According to geostrophic adjustment theory, density 
and pressure lows are associated with cyclonic winds (counterclockwise in the northern hemisphere 
when viewed from above), and density and pressure highs are associated with anti-cyclonic winds 
(clockwise in the northern hemisphere when viewed from above). The plots in Figure 3 indicate that 
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only the dawn-side low for IMF South and the dusk-side low for IMF North are consistent with 
geostrophic force balance. The flows around the low-density cells for IMF East, the dusk-side low for IMF 
South, and the dawn-side low for IMF North are anti-cyclonic. Schoendorf et al. (1996) identify these 
cells as "antibaric" lows. They arise when the drag-driven circulation is strong enough to make 
centrifugal forces non-negligible and capable of producing a low-density cell where a high-density would 
otherwise appear. The results from our simulation at 400 km are consistent with their finding. 


400 km 




(a) 




(c) 



Figure 5. Ion velocity vectors at 400 km (left) and 120 km (right) for (a) IMF South; (b) IMF East; (c) IMF 
North. Vectors are color-coded according to magnitude with black indicating speeds equal or greater 
than 1 km/s for (a) and (b), and 0.5 km/s for (c). 
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Figure 6. Neutral velocity vectors in the rotational frame (top) and neutral density contours (bottom) 
at 400 km. From left to right: IMF East - northern hemisphere; IMF East - southern hemisphere; IMF 
South; IMF North. Vectors are color-coded according to magnitude with black indicating speeds equal 
or greater than 1 km/s. 


Further qualitative validation of the cell structure shown in the density plots of Figure 3 can be obtained 
by comparing them with density maps derived from CFIAMP data by Kwak et al. (2009). Figure 4 
reproduces a plot from their paper. The maps in Figure 4 (Figure 6 in Kwak et al) were generated by 
subtracting the mass density map for zero IMF from mass density maps for each of the four different 
IMF orientations denoted in the figure. By doing this Kwak et al. are able to isolate the effects of IMF 
orientation. Because they are difference maps they are not directly comparable to the plots in Figure 3; 
however, we can still compare relative features like the locations of density highs and lows. Note that 
the Kwak maps were generated over the southern hemisphere and include a case with negative IMF B y . 
While we did not run a negative IMF B y case (the IMF East case had positive IMF B y ) we can use the 
symmetry property that for an IMF with only a B y component, the northern hemisphere ion convection 
pattern for positive IMF B y is the same is the southern hemisphere pattern for negative IMF B y . 
Therefore we expect the leftmost plot of the Kwak figure (IMF B y negative / southern hemisphere) to 
correspond to the leftmost plot of Figure 3 (IMF B y positive / northern hemisphere) and the second to 
left plot of the Kwak figure (IMF B y positive / southern hemisphere) to correspond to the second to left 
plot of Figure 3 (also IMF B y positive / southern hemisphere). While there are differences, many of the 
features seen in the data are reproduced in the ISM simulations. For negative B y (positive B y in the 
northern hemisphere) both data and simulation exhibit enhanced density on the dawn-side, and 
reduced density on the dusk-side (leftmost plots in Figures 3 and 4). This pattern is reversed for positive 
B y (second plot from the left in Figures 3 and 4). For negative B z both data and simulation have a large 
density peak in the cusp region (third plot from the left in Figures 3 and 4). This is the co-located with 
the throat of the convection pattern. The data indicate that the dusk-side low is deeper than the dawn- 
side low, while the ISM simulations have the opposite polarity. It is difficult to compare the case for 
positive IMF B y , since the effects are much weaker for that case (rightmost plots in Figures 3 and 4). 
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Difference Total Mass Density from CHAMP: IMF^0-IMF=0 
B = -4.4 nT B = 4.4 nT B z = -3.4 nT B z = 3.4 nT 



Figure 7. Figure 6 from Kwak et al. (2009). Difference densities at 400 km altitude over the southern 
hemisphere for different IMF orientations. 


Neutral Flows and Density at 120 km 

Figure 5 shows neutral velocity vectors at 120 km altitude for the IMF East, South, and North cases. As in 
the discussion above, we use the southern hemisphere plot for the IMF East case to represent the 
northern hemisphere result from a negative IMF B y simulation. As we observed at 400 km altitude, the 
neutral velocity pattern is characterized by a two-cell convection pattern for the IMF East and South 
cases. For these cases the two-cell pattern is rotated by approximately 4-6 hours to the east 
(counterclockwise in the plots) with respect to the higher-altitude ion convection pattern (see the left 
column of Figure 2). This rotation of the convection pattern is consistent with predictions from the 
analytic model of Mikkelsen and Larsen (1983). It is a consequence of the Coriolis force. We have 
observed in other simulations (not reported on here) that when co-rotation is artificially turned off, the 
velocity pattern aligns with the ion pattern. 

For the IMF North case the flow pattern is more complex, with multiple circulation cells. Closer 
inspection of the neutral velocity pattern for this case indicates that it too is a rotated version of the 
high-altitude ion convection pattern. This is especially evident upon examination of a series of velocity 
vector plots that begin at 400 km altitude and end at 120 km altitude, because the degree of rotation is 
gradual and increases with decreasing altitude. Using this technique we determined that the 120 km 
pattern for the IMF case is rotated approximately 10 hours to the east with respect to the high-altitude 
ion convection pattern. There is also deformation of the circulation cells which is why it is difficult to 
ascertain the degree of rotation in the plot of Figure 5. It is likely that the greater degree of rotation is a 
consequence of the weaker ion driver for this case. 
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Figure 8. Neutral velocity vectors in the rotational frame at 120 km. From left to right: IMF East - 
northern hemisphere; IMF East - southern hemisphere; IMF South; IMF North. Vectors are color-coded 
according to magnitude with black indicating speeds equal or greater than 0.1 km/s. 


Next we examine the neutral density cell structure at 120 km altitude. Figure 6 contains neutral density 
contour plots together with vector plots of the J x B force. Before discussing the rationale for 
comparing / x B vectors with neutral density contours, we will first compare the neutral density cell 
structure to the neutral velocity vectors of Figure 5. For the positive IMF B y case (leftmost plot) we see 
that the density high corresponds to the center of an anti-cyclonic circulation region, and the density 
low centered at morning corresponds to the center of a cyclonic circulation region. Thus these results 
are consistent with geostrophic theory. (Note that the range of the density scale at this altitude is much 
smaller than for the corresponding plots at 400 km). 



Figure 9. JXB vectors (top) and neutral density contours (bottom) at 120 km. From left to right: IMF 
East - northern hemisphere; IMF East - southern hemisphere; IMF South; IMF North. 

Vectors are color-coded according to magnitude with black indicating JXB force 
magnitudes equal or greater than 1X10' 11 dyne/cm 3 . 
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For the surrogate negative IMF B y case (second plot from the left), the density low over the pole 
corresponds to a cyclonic circulation region, but the high is not co-located with the anti-cyclonic 
circulation region. For the IMF South case (third plot from the left) the density lows and highs are 
completely out of phase with their geostrophic circulation centers. And for the IMF North case it is 
difficult to correlate any of the four distinct density cells over the pole with a geostrophic circulation 
center. 

Because of these inconsistencies we looked to an alternative explanation for the 120 km neutral density 
cell structure. Although the neutral flow patterns at 120 km altitude appear to be largely rotational 
convection patterns, as we alluded to in Section 4.2, the ion velocity pattern at 120 km altitude (see 
right column of plots in Figure 2) is largely irrotational with regions of convergent and divergent vector 
directions. This velocity pattern gets impressed upon the neutrals through the ion drag force. As shown 
in Section 2.1 the ion-neutral drag force is proportional to ion mass density, ion-neutral collision 
frequency, and the vector velocity difference between ion and neutral velocities. In steady state it is 
balanced by the J x B force. Thus we look to vector plots of J x B to identify regions where the ion 
drag force is acting to create neutral flow divergences, and thus density enhancements and depressions. 

For all of the cases shown in Figure 6 we see that the regions of convergence correspond to density 
enhancements and regions of divergence correspond to density depressions. Thus the neutral density 
cell structure in regions where the Flail conductivity dominates is organized by the J xB force and 
concomitantly by the rotation of the ion velocity pattern from a predominantly rotational pattern at 
high altitudes to a predominantly irrotational pattern at low altitudes. This is one of the major 
conclusions of this research effort. It suggests that if we know the morphology of the high-altitude 
velocity pattern, and have a reasonable representation of the Flail conductivity, we can predict where 
the high and low density cells will appear. As described in Liu et al. (2004) the prediction of 
pressure/density cells in the thermosphere and their attendant causes could be provide a useful 
framework to order and interpret high-latitude neutral density data. 


4.3. Feedback Mechanisms 

There are a variety of mechanisms by which neutral dynamics can feed back and alter the Poynting flux 
that transmits electromagnetic energy from the magnetosphere to the ionosphere / thermosphere. In 
this section we consider the effects of the drag-driven neutral winds. As momentum is transferred to the 
thermosphere by ion-drag, neutral winds increase and the velocity difference between ions and neutrals 
decreases. The rate at which this occurs depends on collision frequency and neutral mass density. As the 
velocity difference decreases, less current is required from the magnetospheric dynamo (recall from 
Section 2.1 that in steady state the drag force is balanced by the ] xB force). If the dynamo is a 
constant voltage generator, the system responds as it would to a decrease in conductivity, and Poynting 
flux would be reduced. If, on the other hand, the dynamo acts more like a constant current source, the 
system would respond to maintain the velocity difference between ions and neutrals and further 
accelerate both fluids, perhaps resulting in an increase of Poynting flux into the thermosphere. 

To examine these potential effects we can take advantage of the artifice we employed to perform our 
coupled simulations. For each simulation we attained a magnetospheric steady state by first using a 
static neutral atmosphere, i.e. the neutral velocity was set to zero. The thermospheric simulation phase 
was then performed by driving the thermosphere for 24 hours with these fixed magnetosphere 
conditions, allowing the neutral winds to ramp up to equilibrium values. Thus by comparing 
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electrodynamic quantities when neutral velocities are zero (large ion-neutral velocity difference) with 
those after neutral velocities have ramped up (smaller ion-neutral velocity difference) we can assess the 
feedback path discussed above. 

Figure 10 contains plots of electrodynamic quantities for zero neutral velocity (left column) and 
steady-state neutral velocity (right column) for the IMF South case. The top pair of plots contains J x B 
vectors at 300 km, near the F-region conductivity peak. At zero neutral velocity the pattern looks exactly 
like the ion velocity pattern (see the top left plot of Figure 5). This is a consequence of the force balance 
equation of Section 2.1 when neutral velocity is set to zero. For steady-state neutral winds, the 
magnitude of the J x B force is greatly reduced, in agreement with the arguments presented above. 
The bottom pair of plots contains the J x B vectors at 120 km altitude, near the E-region conductivity 
peak. Flere the plots for the two cases are indistinguishable. This is a consequence of the fact that even 
after 24 hours of steady ion-driving, the magnitude of the neutral winds is still a small fraction of the ion 
velocity because of the much larger neutral mass density at this altitude. The middle set of plots shows 
the global effect of the neutral dynamics in terms of the parallel currents that couple the 
magnetosphere with the ionosphere/ thermosphere. These currents are plotted on an altitude surface 
at 500 km. Their close similarity indicates that despite the fact that the F-region is drawing less current, 
the magnitude of this current decrement is so small that it does not alter the global electrodynamics of 
the coupled system. These are controlled at E-region altitudes where the large neutral mass densities 
inhibit neutral flows that would otherwise alter current closure properties. 


5. Conclusions 


We have used a numerical model of the coupled solar wind - magnetosphere - ionosphere - 
thermosphere to elucidate some of the fundamental processes that govern thermospheric response to 
electromagnetic energy inputs. We studied three main subject areas: energy partitioning, horizontal 
neutral density structures, and feedback of neutral dynamics to the coupled system. 

In the area of electromagnetic energy partitioning, we determined that for all of the IMF orientations we 
studied, approximately 90 per cent of the electromagnetic energy goes into Joule heating for 
steady-state conditions. So in a global sense, measurements of Joule heating rates provide a good 
approximation to total electromagnetic energy input; however, locally this is not always the case. For 
example, under southward IMF orientations, local height-integrated Joule heating rates overestimate 
Poynting flux in regions of high conductance, and underestimate Poynting flux in regions of lower 
conductance. This confirms a prediction of Richmond (2010). 

In the area of neutral density structures at F-region altitudes, we confirmed the analysis of Schoendorf 
et al (1996) that attributed the existence of anti-cyclonic low density cells to centrifugal effects in the 
presence of large drag-driven neutral flows. At E-region altitudes, we presented a new theory that 
relates neutral density cell structure to divergences of the ] xB force. At E-region altitudes these 
divergences arise because of Flail conductivity effects which rotate current-density and ion-velocity 
vectors from their non-divergent F-region patterns. 

In the area of feedback effects, we found that increases in neutral wind speeds under steady driving 
conditions do not significantly alter electromagnetic energy flows into the ionosphere / thermosphere. 
While F-region closure currents were significantly reduced, the much larger E-region closure currents 
were not, owing to the much larger neutral densities, and thus much greater coupling times. 
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Figure 10. Electrodynamic quantities for simulations with zero neutral velocity (left column) and 
steady-state neutral velocity (right column), (a) J X B vectors at 300 km altitude; (b) Parallel current 
density contours at 500 km altitude (orange-red indicates current into ionosphere; blue indicates 
current out of the ionosphere); (c)J XB vectors at 120 km altitude 
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7. Publications 


George Siscoe submitted a paper entitled "Aspects of Global Coherence of Magnetospheric Behavior" to 
Journal of Atmospheric and Terrestrial Physics. This paper contains insights gleaned from ISM 
simulations performed for this contract. Contract acknowledgment was included. 

George Siscoe presented a poster at the 2011 GEM/CEDAR meeting titled "Possible explanation of the 
Love-Gannon relationship between Dst and the local-time asymmetry in the low-latitude disturbance 
field". The work he presented was partially funded by this project. 


18. Interactions / Transitions 


Early in this reporting period we consulted with William Burke and Dan Ober of AFRL (Hanscom AFB). 
They provided program direction for the initial phases of this project. 

Keith Siebert gave a presentation at the AFOSR Space Sciences Program Review held in Albquerque, NM 
on 23 June 2011. The title of his presentation was "Electromagnetic Energy Partitioning in the Polar 
Thermosphere". 

Keith Siebert and Jackie Schoendorf met with Dolores Knipp (High Altitude Observatory, NCAR) at the 
ARA Nashua office in July 2011. They discussed Dolores' research interest in intense Poynting flux 
enhancement during IMF B y events (Knipp, 2011). We incorporated her suggestions into our simulation 
analysis. 

Keith Siebert gave an invited presentation at the University of Massachusetts - Lowell Center for 
Atmospheric Research (UMLCAR) on the research described in this report. The title of the presentation 
was "Studies of Electromagnetic Energy Partitioning in the Polar Thermosphere with the Integrated 
Space Weather Prediction Model (ISM)". Interactions with UMLCAR staff provided ideas to achieve 
research objectives for this project. 
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